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Fixed Scale Approach to Equation of State in Lattice QCD 
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A new approach to study the equation of state in finite-temperature QCD is proposed on the 
lattice. Unlike the conventional method in which the temporal lattice size Nt is fixed, the temper- 
ature T is varied by changing Nt at fixed lattice scale. The pressure of the hot QCD plasma is 
calculated by the integration of the trace anomaly with respect to T at fixed lattice scale. This 
"T-integral method" is tested in quenched QCD on isotropic and anisotropic lattices and is shown 
to give reliable results especially at intermediate and low temperatures. 

PACS numbers: 12.38.Gc,12.38.Mh 



I. INTRODUCTION 

The equation of state (EOS) is one of the most funda- 
mental observables to identify different phases in quan- 
tum chromodynamics (QCD) at finite temperature T. 
Also, it is an essential input to describe the space-time 
evolution of the hot QCD matter created in relativistic 
heavy ion collisions So far, the lattice QCD is the 
only systematic method to calculate EOS for wide range 
of T across the region of the hadron-quark phase transi- 
tion (see the recent review, [ij). 

Conventionally, the EOS on the lattice is extracted 
from a method in which T — [Nto)^^ is varied by chang- 
ing the lattice scale a (or equivalently the lattice gauge 
coupling /3 = 6/(7^) with a fixed temporal lattice size Nt- 
Using the thermodynamic relation p — {T/V)\nZ with 
V being the spatial volume and Z beiiig the partition 
function, the pressure p is calculated as j3l| 



T ,^1 dZ T 



dp 



(1) 



Here 5' is the lattice action and (■ • ■) is the thermal 
average with zero temperature contribution subtracted. 
(In multi-parameter cases such as QCD with dynamical 
quarks, should be generahzed to the position vector 
in the coupling parameter space Q.) The initial point 
of integration (3q is chosen in the low temperature phase 
from the condition p{(3o) ~ 0. 

In this conventional method, major part of the com- 
putational cost is devoted to zero temperature simula- 
tions; they are necessary to set the lattice scale and 
to carry out zero-temperature subtractions at the sim- 
ulation points. Furthermore, for the calculation of the 
trace anomaly e — 3p with e being the energy density, the 
non-perturbative beta functions have to be determined 
by zero temperature simulations at the same simulation 
points. In multi-parameter cases, simulations should be 
performed on a line of constant physics (LCP) in the 



coupling parameter space in order to identify the effect 
of temperature to a given physical system. LCP's should 
be determined also at T = 0. These zero temperature 
simulations for a wide range of coupling parameters are 
numerically demanding, in particular, for QCD with dy- 
namical quarks. 

In this paper, we push an alternative approach where 
temperature is varied by Nt with other parameters fixed. 
The fixed scale approach has been applied with the 
derivative method Q in, e.g. Ref. Q. However, determi- 
nation of nonperturbative Karsch coefficients for all the 
simulation points requires quite a work, which will not be 
easy for QCD with dynamical quarks. Here, we propose 
a new method, "the T-integral method": To calculate 
the pressure non-perturbatively, we use 
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which is obtained from the thermodynamic relation valid 
at vanishing chemical potential: 



dT \T^) 
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rp4 
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The initial temperature Tq is chosen such that p{Tq) w 0. 
Calculation of e — 3p requires the beta functions just at 
the simulation point, but no further Karsch coefficients 
are necessary. Since T is restricted to have discrete val- 
ues, we need to make an interpolation of (e — 3p)/T* with 
respect to T. 

Since the coupling parameters are common to all tem- 
peratures, our fixed scale approach with the T-integral 
method has several advantages over the conventional ap- 
proach; (i) T = subtractions can be done by a common 
zero temperature simulation, (ii) the condition to follow 
the LCP is obviously satisfied, and (iii) the lattice scale 
as well as beta functions are required only at the sim- 
ulation point. As a result of these, the computational 
cost needed for T ~ Q simulations is reduced largely. 
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We may even borrow results of existing high precision 
spectrum studies at T = which are pubUc e.g. on the 
International Lattice Data Grid 0. On the other hand, 
when the beta functions arc not available, we need to per- 
form additional T — simulations around the simulation 
points as in the case of the conventional method. 

For a continuum extrapolation, we need to repeat the 
calculation at a couple of lattice spacings. If we adopt 
coupling parameters from T = spectrum studies in 
which the continuum extrapolation has already been per- 
formed, we can get all configurations for T = sub- 
tractions. Furthermore, because the lattice spacings in 
spectrum studies are usually smaller than those used in 
conventional fixed- A^f studies around the critical temper- 
ature Tc, the values of Nt in our approach are much larger 
there than those in conventional studies. For example, 
at a w 0.07 fm, T - 175 MeV is achieved by Nt - 16. 
Therefore, for thermodynamic quantities around T^, we 
can largely reduce the lattice artifacts due to large a 
and/or small Nt over the conventional approach, with- 
out high computational cost for T = calculations. This 
is also a good news for phenomenological applications 
of the EOS, since the temperature achieved in the rela- 
tivistic heavy ion collision at RHIC and LHC will be at 
most up to a few times the critical temperature On 
the other hand, calculation of the trace anomaly around 
and below Tc with large values of Nt may require high 
statistics due to large cancellations by the T — sub- 
traction. We note here that, as T increases, Nt becomes 
small and hence the lattice artifact increases. Therefore, 
our approach is not suitable for studying how the EOS 
approaches the Stephan-Boltzmann value in the high T 
limit. 

Outline of this paper is as follows. After introduc- 
ing our lattice action and the trace anomaly in SectHIl 
we test our T-integral method in SU(3) gauge theory on 
isotropic lattice in Sect lllll and on anisotropic lattice in 
rvl The last section is devoted to summary and conclu- 
sions. 



II. LATTICE ACTION 

We study the SU(3) gauge theory with the standard 
plaquette gauge action on an isotropic and anisotropic 
lattices with the spatial (temporal) lattice size Ng {Nt) 
and lattice spacing (at). The lattice action is given by 



set 


(3 ^ Ns Nt ro/as 


as [fm] 


L[fm] 


a[dg 1 da) 


il 


6.0 1 16 3-10 5.35(13) 


0.093 


1.5 


-0.098172 


12 


6.0 1 24 3-10 5.35(^3) 


0.093 


2.2 


-0.098172 


13 


6.2 1 22 4-13 7.37(3) 


0.068 


1.5 


-0.112127 


a2 


6.1 4 20 8-34 5.140(32) 


0.097 


1.9 


-0.10704 



TABLE L Simulation parameters on Isotropic and anisotropic 
lattices. On Isotropic lattices, we adopt ro/a of [31 , and the 
beta function of Anisotropic ro/as is from [lol |T while the 
beta function is calculated In Sect II VI The lattice scale as 
and lattice size L — NsUs are calculated with ro = 0.5 fm. 
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FIG. 1: Trace anomaly on isotropic lattices as a function 
of physical temperature. The dotted lines are natural cubic 
spline interpolations. Horizontal errors due to the lattice scale 
are smaller than the symbols. 



parameters. The trace anomaly is obtained as 
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rp4 




(6) 



Co \ df3L ie 



(7) 



where ^ = a^/at is the renormalized anisotropy and 
as{dl3/das)^ is the beta function. Note that {d£,o/df3)^ — 
on isotropic lattices. 
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III. EOS ON ISOTROPIC LATTICE 
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where [/^^ (x) is the plaquette in the fiv plane and (3 and 
^0 are the bare lattice gauge coupling and bare anisotropy 



Our simulation parameters are listed in Table [H On 
isotropic lattices, we calculate EOS on three lattices to 
study the volume and lattice spacing dependence. The 
ranges of Nt correspond to T = 210-700 MeV for the sets 

11 and 12, and T = 220-730 MeV for 13, respectively. The 
critical temperature corresponds to Nt ~ 7-8 for il and 

12 and ^ 10 for 13, assuming Tc ~ 290 MeV in quenched 
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FIG. 2; Comparison of the pressure using different interpola- 
tion functions for integration. Statistical errors are estimated 
by the jackknife method [ll|]. 

QCD with the Sommer scale tq = 0.5 fm. The set a2 
will be discussed in the next section. The zero temper- 
ature subtraction is performed with A^t = 16 for il and 
12, and with iV( = 22 for iS. We generate up to a few 
millions configurations using the pseudo-heat-bath algo- 
rithm. Statistical errors are estimated by the jackknife 
analysis. Typically, bin size of a few thousands configu- 
rations are necessary near Tc, while a few hundreds are 
sufficient off the transition region. 

Figure [1] shows (e — 3p)/T'*. Dotted fines in the figure 
are the natural cubic spline interpolations. For compari- 
son, we also reproduce the result of the fixed Nt method 
at A^t = 8 and =32 for which we have rescaled 
the horizontal axis by = 290 MeV according to our 
choice of rg = 0.5 fm. At and below Tc, lattice size de- 
pendence is visible among the sets il {L « 1.5 fm), 12 
(2.2 fm) and the fixed Nt = 8 result (2.7 fm). On the 
other hand, the lattice spacing dependence is negligible 
between il (a w 0.093 fm) and 13 (0.068 fm). At higher T, 
(e — 3p)/T^ on our three lattices show good agreement, 
however they slightly deviate from the fixed Nt one at 
T > 600 MeV, presumably due to the coarseness of our 
lattices at these temperatures. 

The integration of ^ is performed numerically using 
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FIG. 3: EOS on isotropic and anisotropic lattices. Dashed 
horizontal line represents the free gas case in the continuum 
(Stefan-Boltzmann limit). 



the natural cubic spline interpolations shown in Fig. [T] 
For the initial temperature To of the integration, we lin- 
early extrapolate the (e — 3p)/T'* data at a few lowest 
T's because the values of (e — 3p)/T^ at our lowest T 
are not exactly zero. In this study, we commonly take 
To = 150 MeV as the initial temperature which satisfies 
(e — 3p)/T'^ = 0, and estimate the integration from Tq to 
the lowest T by the area of the triangle. 

We estimate the statistical error for the T-integration 
by the jackknife analysis at each T and accumulate the 
contributions from different T's by the error propagation 
because simulations at different A''^ are statistically 
independent. Note that the error for the lattice scale do 
not affect the dimensionless quantity p/T*. Error bars 
shown in the figures represent the statistical errors. 

To estimate the systematic error due to the interpola- 
tion ansatz, we compare the results oip/T^ using cubic 
spline interpolation and those with the trapezoidal rule in 
Fig. [21 (The results on the anisotropic = 4 lattice will 
be discussed later.) We find that the size of systematic 
errors due to the interpolation ansatz are comparable to 
that of the statistical ones. In Fig. [3l we note that the 
natural cubic spline interpolation curve of (e — 3p)/T^ for 
the set 12 shows small negative values at T « 250 MeV 
due to the nearby sharp edge at T « 260 MeV. From a 
comparison with the results of trapezoidal interpolation 
for that range, we find that the effect of this bump on 
the value of p/T^ at T > 250 MeV is 0.032. Although 
a negative pressure is unphysical, because this shift in 
p/T^ is smaller than the statistical errors, we disregard 
the effects of the negative pressure in this paper. To 
avoid arbitrary data handlings, we just adopt the results 
of natural cubic spline interpolations as the central values 
in the foUowings. 

In Fig.[3l we summarize the results of EOS. Results on 
the anisotropic lattice (a2) will be discussed later. The 
normalized energy density e/T^ is calculated by combin- 
ing p/T* and (e-3p)/T4. 
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We find that, except for the vicinity of Tc, the EOS 
is rather insensitive to the variation of lattice size (be- 
tween L 1.5 fm and 2.2 fm) and the lattice spacing 
(between a « 0.093 fm and 0.063 fm). The results of 
EOS agree within 10% for our variation of lattice param- 
eters. This is in part due to the fact that (e — 3p)/T'* 
is not so sensitive to the lattice parameters up to high 
temperatures. Note also that, because (e — 3p)/T'^ = 
in the high temperature limit, the increasingly large lat- 
tice artifacts at large T are naturally suppressed in the 
T-integration. Looking at the T > region closer, we 
note a slight tendency that both p and e decrease as the 
lattice size (lattice spacing) becomes larger (smaller). 

Near and below Tc, we observe a sizable finite size effect 
between L « 1.5 fm (il) and 2.2 fm (12), while the effect 
of the lattice spacing is quite small (il and 13). Therefore, 
for a reliable simulation, we need, at least, L > 2 fm. 

Our results are qualitatively consistent with the pre- 
vious EOS by the fixed Nt method [9J. Quantitatively 
our results are slightly above the fixed Nt results. The 
discrepancy can be in part understood by smaller spatial 
volumes below Tc and the small values of Nt at higher T 
in our method. 



IV. EOS ON ANISOTROPIC LATTICE 

The anisotropic lattice with the temporal lattice finer 
than the spatial one is expected to improve the resolution 
of T without much increasing the computational cost. To 
further test the systematic error due to the resolution of 
T, we perform the study with the T-integral method on 
an anisotropic lattice with the renormalized anisotropy 

The simulation parameters are given as the set a2 in 
Table HI which are the same as those adopted in p^ . 
We vary Nt = 34-8 corresponding to T = 240-1010 
MeV. The zero temperature subtraction is performed 
with Nt = 80. (= 20 x ^). We generate up to a few 
millions configurations. 

We calculate the beta function as{d/3/das)^ b y fi tting 
the ro/tts data with the AUton's ansatz [isj . 

as/ro = R{f3) ■ A (l + Bp{(3) + Cf^P)) , (8) 

where 

/(/3) = i?(/3)/i?(6.10), 5o = 6i = j^^. 

The best fit is achieved by A = 76.0(1.5), B = 0.190(32), 
and C = 0.0204(87) with xVdof = 12.6/6. The beta 
function is now obtained by 

^ (d(3 \ _ I2bl(i l + Bp + Cf^ 
\ da, )^ 6h -boP 1 + 3Bp + 5C/4 ' 
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FIG. 4: Trace anomaly on isotropic and anisotropic lattices 
with similar spatial lattice spacing and lattice size. The dot- 
ted lines are natural cubic spline interpolations. 



and its value at our simulation point is given in Table [D 
For {d^o/dp)^ we adopt the result of [14]. 

In Fig. m we compare the trace anomaly obtained on 
the anisotropic lattice with that on the isotropic lattice 
with similar a, and L (the set 12). We find that the 
results are generally consistent with each other, while, 
due to the cruder resolution in T, the natural cubic spline 
interpolation for the 12 lattice slightly overshoots the data 
on the a2 lattice around the peak at T « 320 MeV. We 
note that the height of the peak on the a2 lattice is similar 
to that of the fine 13 lattice shown in Fig. [T] and the fixed 
Nt = 8 result. Similar to the case of the isotropic lattice, 
our result slightly overshoots the fixed Nt result below 
Tc and at T > 600 MeV. Around the sharp lower edge 
at T « 250 MeV, the undershooting of the natural cubic 
spline interpolation on the 12 lattice, discussed in the 
previous section, is avoided by the finer data points on 
the a2 lattice. 

The results of pressure on the a2 and 12 lattices are 
compared in the lower panel of Fig. [21 While they are 
roughly consistent with each other, we find that the pres- 
sure on the anisotropic lattice is slightly smaller than that 
on the isotropic lattice. This is in part due to the smaller 
values of e — 3p around the peak at T « 320 MeV. We 
also note a slight tendency that the trace anomaly on 
the a2 lattice is systematically lower than that on the 12 
lattice by about la in the high temperature side. 

When we consider the difference between the a2 and 
12 lattices as the systematic error due to the cruder T- 
resolution of the 12 lattice, the systematic error is about 
2-3 times larger than the estimation in the previous sec- 
tion from a comparison of different interpolation ansatze. 
Therefore, for an estimation of the systematic error on 
isotropic lattices, it may be safer to assume a few times 
larger errors than those estimated from a comparison of 
different interpolation ansatze. 

Some part of this difference in the EOS between a2 
and 12 lattices might be explained by the difference of 
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the lattice spacing in the temporal direction, since lat- 
tice artifacts of thermodynamic quantities generally have 
dominant contributions from the temporal lattice spac- 
ing fl^. Indeed this interpretation is supplemented by 
the fact that the difference in the EOS between il and 
13 has a similar tendency as shown in FiglS] We reserve 
further investigation of this possibility for future study. 

V. CONCLUSIONS 

We proposed a fixed scale approach to investigate fi- 
nite temperature QCD on the lattice. To calculate EOS 
non-perturbatively at fixed scale, we introduced the T- 
integral formula ^ , in which p is calculated as an inte- 
gration of the trace anomaly e — 3p with respect to T. 
The fixed scale approach with the T-integral method en- 
ables us to efficiently utilize the results of previous T = 
spectrum studies. At intermediate and low temperatures, 
we can largely reduce the lattice artifacts than the fixed 
Nt approach. On the other hand, our approach is not 
suitable for studying the approach to the high T limit 
because of small Nt there. Also, around and below Tc, 
the large values of Nt may require high statistics due to 
large cancellations by the T = subtraction. 

To test the T-integral method, we performed a series 
of simulations of SU(3) gauge theory on isotropic and 
anisotropic lattices. We found that the method works 
quite well. Our results of EOS for the quenched QCD 
on isotropic and anisotropic lattices are summarized in 
Fig. [31 which are qualitatively consistent with the previ- 
ous results using the conventional fixed N approach. 

Our EOS in the high temperature region turned out 
to be roughly independent of the lattice spacing, lattice 
volume, and the anisotropy. All the results agree within 
about 10% for the range of our lattice parameters. With 
small statistical errors of less than about 2%, we identi- 



fied small systematic shifts under the variations of lattice 
parameters. The smallness of them will be useful for pre- 
cise continuum extrapolations. Around the critical tem- 
perature, wc found that the lattice size should be at least 
larger than about 2 fm. 

The fixed scale approach with the T-integral method 
is applicable to QCD with dynamical quarks too. With 
dynamical quarks, it is more convenient to simulate 
isotropic lattices, because the tuning of anisotropy pa- 
rameters as well as the determination of the factors 
d£,a/d/3 in ([7]) requires quite a few works for full QCD. 
Therefore, it is important to estimate the systematic er- 
ror in EOS due to the limited resolution of T on an 
isotropic lattice. From the test of quenched QCD pre- 
sented in this paper, we found that this systematic error 
is under control and its order of magnitude can be cor- 
rectly estimated by a comparison of different interpola- 
tion ansatze, though it may be safer to introduce a factor 
of about 2-3 to the estimates. 

We are currently investigating EOS in 2+1 flavor QCD 
with non-perturbatively improved Wilson quarks, using 
the configurations by the CP-PACS/JLQCD Collabora- 
tion [l^ , which are public on the ILDG Q . The pseudo- 
critical temperature (~ 175 MeV) is around Nt ~ 16 
on the finest lattice with a w 0.07 fm. We are further 
planning to extend the study to use the 2-1-1 flavor con- 
figurations by the PACS-CS Collaboration generated just 
at the physical values of the light quark masses [13] ■ 
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